# Set working directory to wherever files are saved.
#install.packages('ggeffects')
#library(ggeffects)
#install.packages('ggplot2')
#library(ggplot2)
#install.packages('ggthemes')
#library(ggthemes)

load('introduced_hr_bills.rdata')
load('introduced_senate_bills.rdata')
load('senate_reported_bills.rdata')
load('reported_hr_bills.rdata')
load('engrossed_senate.rdata')
load('engrossed_house.rdata')
load('all_bill_level_agency_scores.rdata')
load('bill_delegation_111.RData') 
load('bill_delegation_110.RData')

#### Figure 2 plots ####

# 110th Congress #

GOP_Bills = colSums(cong110_agency_bills[cong110_agency_bills$Party ==200 & 
                                           cong110_agency_bills$delegates >0 & is.na(cong110_agency_bills$delegates) ==F,
                                         55:ncol(cong110_agency_bills)] !=0)
Dem_Bills = colSums(cong110_agency_bills[cong110_agency_bills$Party ==100 &
                                           cong110_agency_bills$delegates >0 & is.na(cong110_agency_bills$delegates) ==F,
                                         55:ncol(cong110_agency_bills)] !=0)

agency_by_party = cbind(GOP_Bills, Dem_Bills)
agency_by_party = as.data.frame(agency_by_party)
agency_by_party$Agency = NA
agency_by_party$Agency = row.names(agency_by_party)
rownames(agency_by_party) = c()

library(ggplot2)
library(ggthemes)
agency_by_party = gather(agency_by_party, key=Party, value = Bills,
                         c("GOP_Bills", "Dem_Bills"))

agency_by_party$Party[agency_by_party$Party == "GOP_Bills"] = "Republican"
agency_by_party$Party[agency_by_party$Party == "Dem_Bills"] = "Democratic"

agency_by_party$Cabinet = 0
agency_by_party$Cabinet[grep("Department_of", agency_by_party$Agency) ] = 1
agency_by_party$Cabinet[agency_by_party$Agency == "Department_of_the_Air_Force" | 
                          agency_by_party$Agency == "Department_of_the_Army" |
                          agency_by_party$Agency == "Department_of_the_Navy"] = 0

agency_by_party$Proportion = NA
agency_by_party$Proportion[agency_by_party$Party == "Democratic"] = 
  agency_by_party$Bills[agency_by_party$Party== "Democratic"]/sum(agency_by_party$Bills[agency_by_party$Party== "Democratic"])
agency_by_party$Proportion[agency_by_party$Party == "Republican"] = 
  agency_by_party$Bills[agency_by_party$Party== "Republican"]/sum(agency_by_party$Bills[agency_by_party$Party== "Republican"])

agency_by_party$Proportion.c = NA
agency_by_party$Proportion.c[agency_by_party$Party == "Democratic"& agency_by_party$Cabinet == 1] = 
  agency_by_party$Bills[agency_by_party$Party== "Democratic" & agency_by_party$Cabinet == 1]/sum(agency_by_party$Bills[agency_by_party$Party== "Democratic" & agency_by_party$Cabinet == 1])
agency_by_party$Proportion.c[agency_by_party$Party == "Republican"& agency_by_party$Cabinet == 1] = 
  agency_by_party$Bills[agency_by_party$Party== "Republican" & agency_by_party$Cabinet == 1]/sum(agency_by_party$Bills[agency_by_party$Party== "Republican" & agency_by_party$Cabinet == 1])



agencies_party_plot = ggplot(agency_by_party[agency_by_party$Cabinet == 1,], aes(x=Agency, y=Proportion.c, fill=Party)) + #facet_grid(~Party)+
  geom_bar(position=position_dodge() ,stat="identity", color="black") + ylab("Proportion of Bills by Party")+
  theme(axis.text.x = element_text(angle = 90, hjust=1)) + scale_fill_manual(values = c("Blue", "Red"))+ theme_clean() + ggtitle('110th Congress') + scale_y_continuous(limits = c(0, .17))
agencies_party_plot + coord_flip() +theme(legend.position = "bottom")
agencies_party_plot = agencies_party_plot + coord_flip() +theme(legend.position = "bottom") 


# 111th Congress #


GOP_Bills.11 = colSums(cong111_agency_bills[cong111_agency_bills$Party ==200 & 
                                              cong111_agency_bills$delegates >0 & is.na(cong111_agency_bills$delegates) ==F,
                                            55:219] !=0)
Dem_Bills.11 = colSums(cong111_agency_bills[cong111_agency_bills$Party ==100 &
                                              cong111_agency_bills$delegates >0 & is.na(cong111_agency_bills$delegates) ==F,
                                            55:219] !=0)

agency_by_party.11 = cbind(GOP_Bills.11, Dem_Bills.11)
agency_by_party.11 = as.data.frame(agency_by_party.11)
agency_by_party.11$Agency = NA
agency_by_party.11$Agency = row.names(agency_by_party.11)
rownames(agency_by_party.11) = c()


agency_by_party.11 = gather(agency_by_party.11, key=Party, value = Bills,
                            c("GOP_Bills.11", "Dem_Bills.11"))

agency_by_party.11$Party[agency_by_party.11$Party == "GOP_Bills.11"] = "Republican"
agency_by_party.11$Party[agency_by_party.11$Party == "Dem_Bills.11"] = "Democratic"

agency_by_party.11$Cabinet = 0
agency_by_party.11$Cabinet[grep("Department_of", agency_by_party.11$Agency) ] = 1
agency_by_party.11$Cabinet[agency_by_party.11$Agency == "Department_of_the_Air_Force" | 
                             agency_by_party.11$Agency == "Department_of_the_Army" |
                             agency_by_party.11$Agency == "Department_of_the_Navy"] = 0

agency_by_party.11$Proportion = NA
agency_by_party.11$Proportion[agency_by_party.11$Party == "Democratic"] = 
  agency_by_party.11$Bills[agency_by_party.11$Party== "Democratic"]/sum(agency_by_party.11$Bills[agency_by_party.11$Party== "Democratic"])
agency_by_party.11$Proportion[agency_by_party.11$Party == "Republican"] = 
  agency_by_party.11$Bills[agency_by_party.11$Party== "Republican"]/sum(agency_by_party.11$Bills[agency_by_party.11$Party== "Republican"])

agency_by_party.11$Proportion.c = NA
agency_by_party.11$Proportion.c[agency_by_party.11$Party == "Democratic"& agency_by_party.11$Cabinet == 1] = 
  agency_by_party.11$Bills[agency_by_party.11$Party== "Democratic" & agency_by_party.11$Cabinet == 1]/sum(agency_by_party.11$Bills[agency_by_party.11$Party== "Democratic" & agency_by_party.11$Cabinet == 1])
agency_by_party.11$Proportion.c[agency_by_party.11$Party == "Republican"& agency_by_party.11$Cabinet == 1] = 
  agency_by_party.11$Bills[agency_by_party.11$Party== "Republican" & agency_by_party.11$Cabinet == 1]/sum(agency_by_party.11$Bills[agency_by_party.11$Party== "Republican" & agency_by_party.11$Cabinet == 1])



agencies_party_plot.11 = ggplot(agency_by_party.11[agency_by_party.11$Cabinet == 1,], aes(x=Agency, y=Proportion.c, fill=Party)) + #facet_grid(~Party)+
  geom_bar(position=position_dodge() ,stat="identity", color="black") + ylab("Proportion of Bills by Party")+
  theme(axis.text.x = element_text(angle = 90, hjust=1)) + scale_fill_manual(values = c("Blue", "Red"))+ theme_clean()
agencies_party_plot.11 = agencies_party_plot.11 + coord_flip() +theme(legend.position = "bottom") 
agencies_party_plot.11 = agencies_party_plot.11 + scale_y_continuous(limits = c(0, .17))
agencies_party_plot.11 + ggtitle('111th Congress')




#### Table 2 regression models ####

house_sponsors_party_model <- lm(ideol_profile ~ as.factor(Party) + as.factor(Major) + as.factor(cong), data = introduced_hr_bills)

summary(house_sponsors_party_model)

senate_sponsors_party_model <- lm(ideol_profile ~ as.factor(Party) + as.factor(Major) + as.factor(cong), data = introduced_senate_bills)

summary(senate_sponsors_party_model)



#### Creating Figure 3 ####

house_sponsors_model <- lm(ideol_profile ~ DW1 + as.factor(Major) + as.factor(cong), data = introduced_hr_bills)

house_sponsors_model_df <- model.frame(house_sponsors_model)

senate_sponsors_model <- lm(ideol_profile ~ DW1 + as.factor(Major) + as.factor(cong), data = introduced_senate_bills)

senate_sponsors_model_df <- model.frame(senate_sponsors_model)

house_sponsors_model_df$chamber <- 'House'

senate_sponsors_model_df$chamber <- 'Senate'

sponsors_model_df <- rbind(house_sponsors_model_df, senate_sponsors_model_df)

names(sponsors_model_df)[c(3,4)] <- c('Major', 'Congress')

sponsors_ideol_model <- lm(ideol_profile ~ DW1 + Major + Congress, data = sponsors_model_df)

p <- predict_response(sponsors_ideol_model, terms = 'DW1')

sponsors_ideol_plot <- ggplot(p, aes(x, predicted)) + 
  geom_line() + 
  geom_ribbon(aes(ymin = conf.low, ymax = conf.high), alpha = 0.1) +
  labs(x = 'Sponsor DW-NOMINATE', y = 'Bill-Level Agency Ideological Profile') + 
  theme_bw()




#### Table 3 regression models ####

### Figure out what package glm.nb is from (MASS?)

house_cmte_member_model_fullFEs_nb <- glm.nb(delegation_num ~ as.factor(MRef) + as.factor(PooleID) + as.factor(Major), data = introduced_hr_bills, control=glm.control(maxit=100))

summary(house_cmte_member_model_fullFEs_nb)

house_cmte_chair_model_fullFEs_nb <- glm.nb(delegation_num ~ as.factor(ChRef) + as.factor(PooleID) + as.factor(Major), data = introduced_hr_bills, control=glm.control(maxit=50))

summary(house_cmte_chair_model_fullFEs_nb)

senate_cmte_member_model_fullFEs_nb <- glm.nb(delegation_num ~ as.factor(MRef) + as.factor(PooleID) + as.factor(Major), data = introduced_senate_bills, control=glm.control(maxit=100))

summary(senate_cmte_member_model_fullFEs_nb)

senate_cmte_chair_model_fullFEs_nb <- glm.nb(delegation_num ~ as.factor(ChRef) + as.factor(PooleID) + as.factor(Major), data = introduced_senate_bills, control=glm.control(maxit=100))

summary(senate_cmte_chair_model_fullFEs_nb)

#### Creating Figure 4 ####

house_rept <- glm(rept_house ~ delegation_ratio.1 + ImpBill + cong + as.factor(Major) + as.factor(PooleID), data = introduced_hr_bills, family = binomial(link = 'logit'))


senate_rept <- glm(rept_senate ~ delegation_ratio.1 + ImpBill + cong + as.factor(Major) + as.factor(PooleID), data = introduced_senate_bills, family = binomial(link = 'logit'))


house_passed <- glm(PassH ~ delegation_ratio.1 + ImpBill + cong + as.factor(Major) + as.factor(PooleID), data = reported_hr_bills, family = binomial(link = 'logit'))


senate_passed <- glm(PassS ~ delegation_ratio.1 + ImpBill + cong + as.factor(Major) + as.factor(PooleID), data = senate_reported_bills, family = binomial(link = 'logit'))


senate_law_model <- glm(PLaw ~ delegation_ratio.1 + ImpBill + cong + as.factor(Major) + as.factor(PooleID), data = engrossed_senate, family = binomial(link = 'logit'))

house_law_model <- glm(PLaw ~ delegation_ratio.1 + ImpBill + cong + as.factor(Major) + as.factor(PooleID), data = engrossed_house, family = binomial(link = 'logit'))


coef_df <- as.data.frame(matrix(ncol = 5, nrow = 6))

coef_df[1,2] <- house_rept$coefficients[2]
coef_df[2,2] <- senate_rept$coefficients[2]
coef_df[3,2] <- house_passed$coefficients[2]
coef_df[4,2] <- senate_passed$coefficients[2]
coef_df[5,2] <- house_law_model$coefficients[2]
coef_df[6,2] <- senate_law_model$coefficients[2]

coef_df[1,3] <- coef_df[1,2] + 1.96*(summary(house_rept)$coefficients[2,2])
coef_df[1,4] <- coef_df[1,2] - 1.96*(summary(house_rept)$coefficients[2,2])
coef_df[2,3] <- coef_df[2,2] + 1.96*(summary(senate_rept)$coefficients[2,2])
coef_df[2,4] <- coef_df[2,2] - 1.96*(summary(senate_rept)$coefficients[2,2])
coef_df[3,3] <- coef_df[3,2] + 1.96*(summary(house_passed)$coefficients[2,2])
coef_df[3,4] <- coef_df[3,2] - 1.96*(summary(house_passed)$coefficients[2,2])
coef_df[4,3] <- coef_df[4,2] + 1.96*(summary(senate_passed)$coefficients[2,2])
coef_df[4,4] <- coef_df[4,2] - 1.96*(summary(senate_passed)$coefficients[2,2])
coef_df[5,3] <- coef_df[5,2] + 1.96*(summary(house_law_model)$coefficients[2,2])
coef_df[5,4] <- coef_df[5,2] - 1.96*(summary(house_law_model)$coefficients[2,2])
coef_df[6,3] <- coef_df[6,2] + 1.96*(summary(senate_law_model)$coefficients[2,2])
coef_df[6,4] <- coef_df[6,2] - 1.96*(summary(senate_law_model)$coefficients[2,2])

coef_df[1:2,1] <- 'Reported from Committee'
coef_df[3:4,1] <- 'Passed Chamber'
coef_df[5:6,1] <- 'Became Law'

coef_df[c(1,3,5),5] <- 'House'
coef_df[c(2,4,6),5] <- 'Senate'

names(coef_df)[c(1:5)] <- c('model', 'deleg_ratio_coef', 'upper_CI', 'lower_CI', 'Chamber')


coef_df$model <- c('Reported from Committee', 'Reported from Committee', 'Passed Chamber', 'Passed Chamber', 'Became Law', 'Became Law')

coef_df$model <- factor(coef_df$model, levels=c("Reported from Committee", "Passed Chamber", "Became Law"))

p <- ggplot(coef_df) + geom_pointrange(position = position_dodge(0.5), aes(model , deleg_ratio_coef, ymin = lower_CI, ymax = upper_CI, color = Chamber)) + geom_hline(yintercept = 0) + theme_bw()

p <- p + xlab('Stage') + ylab('Delegation Ratio Coefficient')

p


#### Table 4 regression models ####

hr_bills_model <- betareg(delegation_ratio.1 ~ as.factor(cong) + house_median_dist + as.factor(Major) + as.factor(version) + dim1 + dim2, data = subset(all_bill_level_agency_scores, all_bill_level_agency_scores$bill_type == 'hr' & all_bill_level_agency_scores$delegation_num > 0))

summary(hr_bills_model, type = 'deviance')

hr_bills_passed_house_model <- betareg(delegation_ratio.1 ~ as.factor(cong) + house_median_dist + as.factor(Major) + as.factor(version) + dim1 + dim2, data = subset(all_bill_level_agency_scores, all_bill_level_agency_scores$bill_type == 'hr' & all_bill_level_agency_scores$delegation_num > 0 & all_bill_level_agency_scores$PassH == 1))

summary(hr_bills_passed_house_model, type = 'deviance')

senate_bills_model <- betareg(delegation_ratio.1 ~ as.factor(cong) + senate_median_dist + as.factor(Major) + as.factor(version) + dim1 + dim2, data = subset(all_bill_level_agency_scores, all_bill_level_agency_scores$bill_type == 's' & all_bill_level_agency_scores$delegation_num > 0))

summary(senate_bills_model, type = 'deviance')

senate_bills_passed_senate_model <- betareg(delegation_ratio.1 ~ as.factor(cong) + senate_median_dist + as.factor(Major) + as.factor(version) + dim1 + dim2, data = subset(all_bill_level_agency_scores, all_bill_level_agency_scores$bill_type == 's' & all_bill_level_agency_scores$delegation_num > 0 & all_bill_level_agency_scores$PassS == 1))

summary(senate_bills_passed_senate_model, type = 'deviance')

law_model <- betareg(delegation_ratio.1 ~ as.factor(cong) + cong_median_dist + as.factor(Major) + dim1 + dim2, data = subset(all_bill_level_agency_scores, all_bill_level_agency_scores$version == 'enr' & all_bill_level_agency_scores$delegation_num > 0))

summary(law_model, type = 'deviance')

mayhew_model <- betareg(delegation_ratio.1 ~ as.factor(cong) + cong_median_dist + dim1 + dim2, data = subset(all_bill_level_agency_scores, all_bill_level_agency_scores$version == 'enr' & all_bill_level_agency_scores$delegation_num > 0 & all_bill_level_agency_scores$mayhew == 1))

summary(mayhew_model, type = 'deviance')
